rm(list = ls())

library(tidyverse)
library(lfe)

## Load MP data

dat <- readRDS('data/mp_data.rds')  %>%
  filter(vote_abstain_no_show == 0) %>%
  filter(elecper != 8) 

experience_covars <- c('nsdap_member', 
                       'syn_in_gem_bin', 
                       'rel_cath',
                       'veteran_ww1',
                       'veteran_ww2', 
                       'capture_ww2',
                       'soviet_capture',
                       'exile_repressed_combined', 
                       'resistance_member')

experience_covars_labels <- c('NSDAP member', 
                              'Jewish presence (synagogue)',
                              'Catholic',
                              'WW1 Veteran', 
                              'WW2 Veteran',
                              'WW2 POW',
                              'Soviet POW',
                              'Repressed / Exile', 
                              'Resistance')


covars <- c('year_birth', 'gender', 'mandate', 'dualcand', 'closeness_district_categ')

##

m1 <-  paste('vote_yes', '~', paste(experience_covars[1:3], collapse = "+"),'', ' | elecper  | 0 | ags_cluster')
m2 <-  paste('vote_yes', '~', paste(experience_covars[1:3], collapse = "+"),'', ' | elecper + party_fe | 0 | ags_cluster')
m3 <-  paste('vote_yes', '~', paste(experience_covars[1:3], collapse = "+"),'', ' | elecper + party_fe + state_id | 0 | ags_cluster')
m4 <-  paste('vote_yes', '~', paste(experience_covars, collapse = "+"),'+', paste(covars, collapse = "+"), ' | elecper  | 0 | ags_cluster')
m5 <-  paste('vote_yes', '~', paste(experience_covars, collapse = "+"),'+', paste(covars, collapse = "+"), ' | elecper + party_fe | 0 | ags_cluster')
m6 <-  paste('vote_yes', '~', paste(experience_covars, collapse = "+"),'+', paste(covars, collapse = "+"), ' | elecper + party_fe + state_id | 0 | ags_cluster')

mlist <- list(m1, m2, m3, m4, m5, m6)

## Estimate 

mlist_res <- lapply(mlist, function(m){
  
  m_est <- felm(m %>% as.formula(),
                data = dat)
  
  return(m_est)
  
})

## To Table 

stargazer::stargazer(mlist_res, 
                     style  = 'ajps',
                     keep = c('rel_cath', 'nsdap_member', 'syn_in_gem_bin'),
                     covariate.labels = c('NSDAP member', 'Jewish presence', 'Catholic'))

## Plot for main spec

m_main <- mlist_res[[6]] %>%
  broom::tidy(conf.int = T) %>%
  filter(term %in% experience_covars) %>%
  mutate(label = experience_covars_labels,
         label = fct_reorder(label, estimate)) 

m_main

p1 <- ggplot(m_main, aes(x = label, y = estimate, ymin = conf.low, ymax = conf.high)) + 
  geom_errorbar(width = 0) + 
  geom_point(shape = 21, fill = 'white') + 
  theme_bw() + 
  coord_flip() + 
  geom_hline(yintercept = 0, linetype = 'dotted') + 
  labs(x = '',
       y = 'Effect estimate')

p1











